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We describe a large scale simulation effort using OPAL (Object Oriented Parallel Accelerator 
Library), that leads to a better quantitative understanding of the existing PSI high power proton 
cyclotron facility. The 1.3 MW of beam power on target poses stringent constraints on the controlled 
and uncontrolled beam losses. We present initial conditions for the Ring simulation, obtained 
from the new time-structure measurement and the many profile monitors of the 72 MeV transfer 
line. A trim coil model is developed, needed to avoid the dangerous Vr — 2vz resonance. By 
properly selecting the injection position and angle (eccentric injection), the flattop voltage and 
phase, very good agreement between simulations and measurements at the radial probe RRE4 is 
obtained. We report on 3 ... 4 orders of magnitude in dynamic range when comparing simulations 
with measurements. The relation between beam intensity, rms beam size, and accelerating voltage 
is studied and compared with measurement. The demonstrated capabilities are mandatory in the 
design and operation of the next generation high power proton drivers. In an outlook we discuss 
our future plans to include more physics into the model, which eventually leads to an even larger 
dynamic range in the simulation. 

PACS numbers: 29.20.dg;29.27.Bd;41.85.Ew 

I. INTRODUCTION 

PSI operates a cyclotron based high intensity proton accelerator routinely at an average beam 
power of 1.3 MW. With this power the facility is at the worldwide forefront of high intensity proton 
accelerators. An upgrade program is under way to ensure high operational reliability and push the 
intensity to even higher levels. The beam current is limited in practice by losses at extraction and 
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the resulting activation of accelerator components. Further intensity upgrades and new projects 
aiming at an even higher average beam power, are only possible if the relative losses can be lowered 
in proportion, thus keeping absolute losses at a constant level. 

Maintaining beam losses at levels allowing hands-on maintenance is a primary challenge in any 
high power proton machine design and operation. For a 1.3 MW beam in the PSI Ring cyclotron 
this corresponds to a transmission of 99.97% taking controlled and uncontrolled losses into account. 
In a 10 MW class machine we require the losses to be on the same level which is a challenging 
task and is asking for precise beam dynamics calculation. In consequence, predicting beam halo 
at these levels is a great challenge and will be addressed in this paper. 

High power hadron drivers have being used in many disciplines of science and, a growing interest 
in cyclotron technology for high power hadron drivers has be shown very recently. Two very recent 
papers demonstrate this fact: 1) The search for CP violation in the Neutrino sector [1] calls 
ultimately for three machines in the megawatt range at an energy of 800 MeV. 2) in [2], a white 
paper on Accelerator and Target Technology for Accelerator Driven Transmutation and Energy 
Production the cyclotron technology is advertised, quote: On the whole, the development status of 
accelerators is well advanced, and beam powers of up to 10 MW for cyclotrons and 100 MW for 
linacs now appear to be feasible 

This report will briefly introduce OPAL, a tool for precise beam dynamics simulations including 
3D space charge. One of OPAL's "flavors" (OPAL-cycl) is dedicated to high power cyclotron 
modeling and is explained in greater detail. We then explain how to obtain initial conditions for 
our PSI Ring cyclotron which still delivers the world record in beam power of 1.3 MW in continuous 
wave (cw) operation. Several crucial steps are explained, necessary to be able to predict tails at 
the level of 3cr . . . 4cr in the PSI Ring cyclotron. We compare our results at the extraction with 
measurements, obtained with a 1.18 MW cw production beam. Based on measurement data, 
we develop a simple linear model to predict beam sizes of the extracted beam as a function of 
intensity and confirm the model with simulations. A conclusion and discussions to include more 
physics into the model, which eventually leads to a even larger dynamic range in the simulation, 
closes the paper. 
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II. BASIC EQUATIONS AND PHYSICAL MODEL 
A. A BRIEF LOOK AT OPAL 

OPAL (Object Oriented Parallel Accelerator Library) is a tool for charged-particle optic cal- 
culations in large accelerator structures and beam lines including 3D space charge. OPAL is 
built from first principles as a parallel application, it admits simulations of any scale: on the 
laptop and up to the largest High Performance Computing (HPC) clusters available today. Sim- 
ulations, in particular HPC simulations, form the third pillar of science, complementing theory 
and experiment. OPAL includes various beam line element descriptions and methods for single 
particle optics, namely maps up to arbitrary order, symplectic integration schemes and lastly time 
integration [3]. OPAL is based on IPPL (Independent Parallel Particle Layer) [4] which adds 
parallel capabilities. Main functions inherited from IPPL are: structured rectangular grids, fields 
and parallel FFT and particles with the respective interpolation operators. Recently we added a 
powerful iterative solver to OPAL taking into account complicated boundary conditions [5]. More 
details on cyclotron modeling which are direct relevant to this article can be found in [6] . Several 
flavors of OPAL are available. For details we refer to the User Manual [3]. In this paper we use 
OPAL-t for the tracking of 72 MeV beam line, connecting two cyclotrons, the Injector 2 and the 
Ring Cyclotron. The other OPAL flavor - OPAL-cycl - is designed specially for cyclotron beam 
dynamics and, is explained in the next section. 



B. THE BEAM DYNAMICS MODEL OF OPAL-cycl 

In the cyclotrons and beam lines under consideration, the collisions between beam particles can 
be neglected because the typical bunch densities are low. In time domain, the general equations 
of motion of a charged particle in electromagnetic fields can be expressed by 

^=g(c^xB + E), (1) 

where mo, 9, 7 are rest mass, charge and the relativistic factor. With p — m^c^(3 we denote the 
momentum of a particle, c is the speed of light, and (3 ~ (Pxy f^y, f^z) is the normalized velocity 
vector. In general the time (t) and position (x) dependent electric and magnetic vector fields are 



written in abbreviated form as B and E. 

If p is normalized by toqC, Eq. (1) can be written in Cartesian coordinates as 

at ttiqc jiTiQ 

^ - -^Ey + ^ip^B^-pM, (2) 

dt TUqC JITT-q 

-TT = E^-l [p^By - PyB^}. 

at tuqc ^rriQ 

The evolution of the beam's distribution function /(x, cf3, t) can be expressed by a coUisionless 
Vlasov equation: 

^^^dtf + cf3- V,/ + q{E + c/3 X B) • V^^/ = 0, (3) 
where E and B include both external applied fields, and space charge fields 

E = Ecxt + Esc, 

B = Bext+Bsc. (4) 

In order to model a cyclotron, the external electromagnetic fields are given by measurements or 
by numerical calculations. 

The space charge fields can be obtained by a quasi-static approximation. In this approach, 
the relative motion of the particles is non-relativistic in the beam rest frame, so the self-induced 
magnetic field is practically absent and the electric field can be computed by solving Poisson's 
equation 

V^</>(x) = -4^, (5) 
So 

where cj) and p are the electrostatic potential and the spatial charge density in the beam rest frame. 
The electric field can then be calculated by 

E = -V0, (6) 

and back transformed to yield both the electric and the magnetic fields, in the lab frame, required 
in Eq. (4) by means of a Lorentz transformation. Because of the large vertical gap in our cyclotron, 
the contribution of image charges and currents are minor effects compared to space charges [7], 
and hence it is a good approximation to use open boundary conditions. 
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FIG. 1: (Color) Schematic plot of the top view of the three coordinate frames. The red curve is the orbit of 
the bunch center, the blue area represents the bunch shape, and the gray area is the hill region of magnetic 
field. 

The combination of Eq. (3) and Eq. (5) constitutes the Vlasov-Poisson system. In the following, 
the method of how to solve these equations in cyclotrons using PIC methods is described in detail. 

Considering that particles propagate spirally outwards in cyclotrons, and the longitudinal orien- 
tation changes continuously, three right-handed Cartesian coordinate systems are defined, as shown 
in Fig. 1. The first coordinate system is the fixed laboratory frame Siab, in which the external field 
data is stored and the particles are tracked. 

Its origin is the center of the cyclotron and its X ~Y plane is the median plane and the positive 
direction of Z axis points vertically upwards. 

The second coordinate system is the local instantaneous frame Siocai, which is a temporal aux- 
iliary frame for the space charge solver. Its origin O' is the mass center of the beam and the 
orientation of the Y' axis is coincident with the average longitudinal direction and the positive 
orientation of the Z' axis points vertically upwards. 

The third coordinate system is the beam rest frame Sbcam, which is co-moving with the centroid 
of the beam. It has the same orientation and origin as Siocai, but the length in longitudinal 
direction is scaled by I/7 due to relativistic effects. 

At each time step, the frames Siocai and Sbcam are redefined according to the current 6D phase 
space distribution, and all particles are transformed from Slab to Siocai , then a Lorentz transfor- 



mation is performed to transform all particles to Sbcam- The Poisson equation is then solved in 
the frame Sbcam- In ^ 3D Cartesian frame, the solution of the Poisson equation at point [x^y,z) 
can be expressed by 

'/'(a;, y, z) = j G{x, x', y, y' , z, z')p{x', y' , z')dx'dy'dz' , (7) 

with G the 3D Green function 

G(x, x' ,y,y' z, z') = — , (8) 

^{x~x')^ + {y^y'r + (z-z')^' 

assuming open boundary conditions. Details of the space charge field calculation can be found in 

[8]. 

The model of the external magnetic field is based on mid-plane field measurements with excited 
trim coils. In consequence we have a vertical field, Bz, measured on the median plane (z = 0) as a 
function of azimuthal position (0). Since the magnetic field outside the median plane is required to 
compute trajectories with z ^ 0, the field needs to be expanded in the Z direction. According to 
the approach given by Gordon and Taivassalo [9] , by using a magnetic potential and the measured 
Bz on the median plane at the point {r,9,z) in cylindrical polar coordinates, the 3rd order field 
can be written as 

Br{r,e,z) = z^-lz^Cr, 
or 6 

Be{r,e,z) = '-^-^'-Ce, (9) 
r aO 6 r 

B,{r,e,z) = Bz-\z^Gz, 

where B^ = Bz{r,9,0) and 

C = ^'-g^ I 1 9^Bz 1 dBz ^ 1 d^Bz ^ 1 d^Bz 

^ dr^ r dr"^ dr drd6^ 89'^ ' 

I d^Bz d^Bz 1 d^Bz 

^ ~ r drdO ^ dr^de ^ 86^ ' ^ ' 

_ IdBz d^Bz 1 d^Bz 
Gz — 



r dr dr"^ r"^ dO"^ 
All the partial differential coefficients are computed on the median plane data by interpolation, 
using Lagrange's 5-point formula. 

Finally both the external fields and space charge fields are used to track particles for one time 
step using a Ath order Runge-Kutta (RK) integrator, in which the fields are evaluated for four 
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O locations of 3 time-structure probes 




FIG. 2: (Color) The PSI Injector 2 cyclotron and beam transfer line to the PSI Ring cyclotron. 

times in each time step. Space charge fields are assumed to be constant during one time step, 
because their variation is typically much slower than that of external fields. More details and 
unique features can be found in [6]. 

III. OBTAINING INITIAL CONDITIONS FOR THE RING CYCLOTRON 

At the extraction region of the Injector 2 we only have a very limited number of measurement 
data, however in the injecting beam line connecting the two cyclotrons we have 14 vertical and 
17 horizontal beam profile monitors available for high intensity operation. Three time-structure 
measurements, one at the last turn of the Injector 2, one 27 meters downstream and one at the 
first turn of the PSI Ring cyclotron give important information on the longitudinal beam size [10]. 
In an overview (Fig 2) the starting point of the simulations, and some of the diagnostics are shown. 
We note that from a beam dynamics point of view, the particles travel in the order of 4 km, from 
the marked start of the simulation to the RRE4, the probe covering the last 9 turns of the PSI 
Ring cyclotron. 

We start the OPAL-t simulations from the middle of the last valley before extraction from the 
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Injector 2 and perform a full 3D simulation until the magnetic injection channel (MIC) of the 
PSI Ring cyclotron. In Tab. I the initial values for the simulation of important beam parameters 
are shown. At the MIC we resample the distribution and switch to OPAL-cycl for the PSI Ring 
cyclotron simulation. The new distribution is sampled using the moments obtained from the 
transfer line simulation (at MIC). 

TABLE I: Initial conditions of the 72 MeV transfer line for a 2 mA cw beam. The emittances are non- 
normalized. 

Distribution Ex (mm-mrad) ey (mm-mrad) Xj-ms (mm) j/r-ms (mm) Irms (mm) S (%) {xx') {yy') {x5) 
3D Gaussian 2.22 0.43 2.9 0.5 6.2 0.06 -0.14 0.07 -0.92 

Figure 3 shows the comparison of the 2(t beam width between OPAL-t and the measurements. 
The model predicts very well the evolution of the envelope from the beginning to the end of the 
transfer line. 
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FIG. 3: Envelope of the beam in the 72 MeV transfer line for a 2 mA beam. The error bars showing 
different measurements at the same intensity level 
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Figure 4 shows the comparison of the predicted bunch length by the model and the measurements 
using the time-structure probes. The large error bar at the 

injection to the Ring is because of the large background. The longitudinal initial conditions, in 
the center of the valley between sector magnet 3 and 4 (see Fig: 2), are derived from the time 
structure measurement, which is ^ 4 m upstream. In the center of the valley, the major axis of 
the bunch ellipse is along the longitudinal direction, and hence we obtain the transformed initial 
condition by a simple rotation in the longitudinal and radial plane. 




FIG. 4: Bunch length of the beam in the 72 MeV transfer line for a 2 mA beam (in Fig 2 the positions of 
the time-structure measurements are shown). 



IV. TOWARDS REALISTIC HIGH POWER CYCLOTRON SIMULATIONS 

The beam losses during the operation of the cyclotron usually limits the intensity that can be 
extracted. The PSI 590 MeV Ring routinely delivers 2.2 mA of cw beam, having a very low 
integrated loss rate, of the order of 0.02%. This thight margin avoids excessive activation of 
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accelerator components and hence keeps the radiation dose imposed on the personnel involved in 
maintenance at acceptable levels. Furthermore, about 90% of the losses are located at injection 
and extraction. 

Therefore, the understanding of the beam dynamics and, the knowledge of the detailed beam 
distribution especially at the extraction region, is one of the key points to be addressed especially 
if power levels increase in future projects [2]. 

Several important effects which need to be carefully modeled to keep extraction losses in the 
order of 0.02% are: 

• the turn separation at the position of the extraction septum must be made as large as possible, 

• the radial beam size at the extraction region must be smaller than the turn separation, 

• the halo, especially at the extraction, has to be minimized, 

• in case of the PSI Ring cyclotron, a long "pencil" beam is used and hence the linear space 
charge effects must be effectively compensated to avoid the formation of a S-shaped beam 
which apparently increases the effective radial beam size. 

We now discuss these issues related to the PSI Ring cyclotron which however can be considered 
universal for high power cyclotrons, and hence are certainly important for future high intensity 
related projects [2]. 

A. The Flattop Phase 

Although a compact beam is observed at the extraction of the Injector 2 cyclotron, the bunch 
length increases from about a = 6 mm to about a = 31 mm at injection into the Ring after passing 
through the almost 60 m long (72 MeV) transfer line. For such a long "pencil" beam, a flattop 
cavity is needed to compensate the energy difference from the main cavity and avoid the formation 
of the S-shape beam caused by space charge effects. 

When there is no space charge effect, the ideal flattop makes the total energy gain of any 
particle almost the same independent of the RF phase. Considering a high current beam, the 
flattop phase must be shifted such that the tail particles gain more energy than the head particles. 
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This compensates exactly the hnear part of the space charge force. Therefore the phase of the 
flattop is adjusted intensity-dependent and, there exists an optimum flattop phase for a given 
intensity. 

For our simulation we use 11.5% of the sum of the main cavity voltages as the flat-top cavity 
voltage (as set in the control room) and adjust the phase to obtain the same phase difference 
between main and flat-top cavity as set in the control room. 



B. The Effect of the Trim Coil TCI 5 



A small manufacturing error produces a slight deviation in the average field profile and a cor- 
responding shift in the tunes Vr and v^. This requires a strong excitation of the trim coil TC15. 
Without this correction, the coupling resonance — "ivz would be crossed four times at energies 
of 490, 525, 535 and 585 MeV, respectively. A large horizontal oscillation would be transformed 
into a large vertical one at the coupling resonance which can lead to large vertical beam losses. An 
analytic model was developed which mimics the field due to real trim coil characteristics [11]. It 
is described by Eq. (11). 

A2 , A3 



AS = -B„ 



(11) 



where R is the radius, Ba is the maximum magnetic field, the constants Ai = —1.08, A2 = 1.08, 
^3 = 1.80, ai = 0.005, bi = -21.72, 02 = -0.033, 62 = 145.60 for our case. It provides an 
additional magnet field and field gradient in the radial direction as shown in Fig. 5. 
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FIG. 5: The field (solid line) and field gradient (dashed line) of TC15 in the PSI Ring cyclotron. 



The trim coil gives a maximum magnetic field of 14 Gauss, and furthermore has a long tail 
towards smaller radii in order to make the integrated strength of the trim coil over the radius to 
zero. The radial and vertical tune shifts caused by TC15 are given by, 



_J3_dB Q Q-,^ 



(12) 



-0.028 



where R is the orbit radius, B is the hill field, ^ is the average field gradient in radial direction. 
Careful beam dynamics studies have shown the meaningfulness of such detailed modeling in order 
to obtain a complete and precise pictures of the beam dynamics in the PSI Ring cyclotron. The 
modified tune diagram by TC15 is shown in figure. 6. 
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FIG. 6: Tune diagram with and without TC15. 



Without TC15, in simulations and in the operation of the Ring, we observe severe vertical beam 
losses and can not obtain the required extraction efficiency. 



C. The Injection Position and Angle 

The PSI Ring cyclotron has a single turn extraction, hence a large radial turn separation between 
the last two turns is required. The turn separation for a centered beam is defined as 

dn + 1 E 1 + k ^ ' 

where k is the field index. For the PSI Ring cyclotron this gives about 6.0 mm (Fig. 7 upper part) 
at the extraction region, which is not enough for high current operation and would result in large 
losses. 

To increase the turn separation, a non-centered injection into the PSI Ring cyclotron is used. 
Since Vr ~ 1-7 at extraction, adjusting the injection position and angle, results in the betatron 
amplitude being almost equal to the increase in radius per turn. The formation of the turn pattern 
under this condition, for the last nine turns, is shown in Fig. 7 (lower part). 
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FIG. 7: Schematic representation of the turn pattern in the PSI Ring cyclotron. Ic, (l>c and Ir represent 
the intensity distribution of a centered beam with 6 mm turn separation, the betatron oscillation phase of 
an eccentric beam and the intensity distribution of a real beam with eccentric injection. 

This is a special turn pattern because the last turn is well separated from the overlapping second, 
third and fourth last turns. In this case, the turn separation at the extraction turn is as large as 
16 mm, hence it allows the extraction of a high intensity beam with very low losses. 

D. Comparing the Radial Intensity Profile at Extraction with Measurements 

Up to now we have described the most important steps in setting up a precise beam dynamics 
simulation of the PSI Ring cyclotron. We now compare simulations with measurements from a 
radial probe (RRE4) covering the last 9 turns of the PSI Ring cyclotron. The probe is located 30 
cm upstream from the 50 /im thick electrostatic extraction septum and, hence gives a very good 
picture of the beam distribution at the septum. 

This probe is able to measure at the full intensity of the 1.3 MW cw beam. In order to compare 
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the simulations with measurements, not only a radial probe is implemented in OPAL-CYCL but 
also all other parameters, described in the previous sections, can be entered into the simulation. 

The flattop phase and the injection position and angle are optimized to get the largest turn 
separation and smallest beam size at the extraction region, in both the simulation and operation 
of the PSI Ring cyclotron. 

The effect of the trim coil TC15 on the turn pattern is shown in Fig. 8, black denotes the 
measurement and the colors distinguishing simulations with and without TC15 . For a fixed 
energy the shift is given by = For turn 180 this shift is: AR \max~ 3mm, hence the 

center of turn 180 moves to the exact position of the measurement when considering the effect of 
the TC15. 

^ ' U82 ' ^185 ' ' '■ 




r (mm) 

FIG. 8: Radial beam profile with indicated turn numbers at extraction for a 2 mA beam, for a parabolic 
initial distribution at MIC 

In the PSI cyclotron facility, the beam is heavily collimated during the early stage of acceleration 
in the Injector 2 and in the beam transfer line to the PSI Ring cyclotron. As a consequence, 
the beam profiles do not follow a Gaussian distribution. In Fig. 9 we show again the intensity 
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pattern at the last 9 turns but for different starting distributions at MIC. Using the properties of a 
binomial distribution [12], we can vary a single parameter m from to cx) and cover a wide range 
of distributions: from KV to Gaussian. We find that a parabolic distribution (m — 2) matches 
best the measurement at the crucial point, the septum. On the extreme side, as expected, the 
Gaussian distribution (truncated at Sti) with its tails would fill up the intensity dip and hence 
would increase the losses at the septum. This indicated that there are indeed sharp edges in the 
real distribution, which still differs from the assumed idealized distribution, as suggested by the 
remaining differences between simulation and measurement (Fig 9, e.g at i? = 4434). 




r (mm) 



FIG. 9: The comparison of different initial beam distributions at the extraction probe, for a 2 mA beam. 

Nevertheless, this remarkable agreement is obtained after « 4km of tracking the millions of macro 
particles through external- and self fields. This is only possible because of the parallel nature of 
OPAL which allows such simulations to be carried out on large high performance computing 
clusters. The statistical error of the measurement is indicated at i? = 4423 mm in Fig. 8 and 
9. At other radii the error bars are significantly smaller and hence not shown in the figures. The 
statistical errors of the simulation are smaller than those of the measurement errors due to the 
large number of simulation particles. 
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E. SCALING LAW OF BEAM SIZE WITH RESPECT TO CURRENT 

The energy spread AEsc{linear) — eAUsc{linear) caused by the Hnear longitudinal space charge 
force after n revolutions is given by Joho [13] using the sector model: 



where gic ~ 1.4 (form factor), Zq — 1/Att£qc = 30i7 and (/) is the average current, /S.<f) is the phase 
width, n the turn number and /3y — Vj j where is the final velocity of the beam. 

The linear energy spread can be compensated with a tilted flattop voltage. This reduces the 
energy of leading particles and increases the energy of trailing particles. There remains a non 
linear part of the energy spread AEsdnonlinear) , which can not be compensated. Let's define 



and note that /„ depends strongly on the beam distribution and is in our case an open parameter 
in the range of 0.1 . . . 0.5. 

According to Eq. 13, the space charge induced energy spread leads to a radial spread (Ai?sc) 
that results in an increase in beam size. In Fig. 10 we compare beam sizes at the extraction for 
beam currents from 10 fiA to 2.2 mA with simulations. 



AC/, 



SC 




(14) 



AEsc{nonlinear) = fn AEsdHnear), 



(15) 
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FIG. 10: Radial beam size at extraction vs. beam current. 



Even though the measurements where done over a time span of 4 years with very different 
machine configurations we obtain a good agreement between the simulations (theory) and the 
measurements. Hence we can predict very well the extracted beam size as function of the intensity. 
This simple model and the precise simulations shown in the previous paragraph constitutes a 
benchmarked model for the prediction of the most delicate parameters in high intensity cyclotrons. 

The relation between the average energy gain dE/dn and the turn number n is: 



n'^=Ef~E,. (16) 
dn ' ^ ' 

where Ef is the final energy and Ei is the initial energy. For single turn extraction the loss on the 
septum is limited by the ratio 



AEsc (nonlinear) 



dE 
dn 
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leading to the condition 

dE 

AEscinonlinear) < fin—;—- (17) 
an 

We obtain empirically, for a centered beam a value for /i„ which is fin ~ l/^j where as for an 
eccentric beam, the turn separation is enhanced, and hence fin can be as high as fin ~ 1- 

Putting (14), (15), (16) and (17) together, we get for the current limit from longitudinal space 
charge forces 

Ki max — J. ry So' V-*-"/ 

where Uf — Ef/e and Ui = Ei/e. Since the turn number n is inverse proportional to the cavity 
voltage Vcav, we see the big advantage of a large cavity voltage 



{I)max ^ ^3 ^ ^cav (19) 

as experimentally demonstrated by the historical current development in the PSI Ring cyclotron 
shown in Fig 11. 
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FIG. 11: Production beam current versus turn number over a span of 32 years of the PSI Ring cyclotron 
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This relation was first predicted by Joho [13] and confirmed during almost 36 years of operation 
of the PSI Ring cyclotron. 

The prediction of the current limit in the PSI Ring cyclotron gives with (18) {I)max = 3mA, 
using the following parameters: = 72 MV, Uf = 590 MV, /?/ = 0.8, n = 188 turns {dE/dn = 3 
MeV), A(/) = 12 and estimating /„ = 1/4 and /i„ = 1 (eccentric injection). 

This is remarkably close to the present current limit of 2.3mA (2010), given the crude assump- 
tions: 

• no turn structure inside the charge sheet (see [13] figure 3), 

• non relativistic approximation, 

• no radial boundary condition and 

• uncertainties in /„ and /l(„. 

We note that for the PSI Injector 2 (18) is not applicable due to phase mixing [14] and hence would 
give a pessimistic value for the current limit. 

V. CONCLUSIONS AND DISCUSSIONS 

In this paper, we present novel precise simulations for the beam dynamics in high intensity 
cyclotrons. For the first time we are able to obtain a realistic and detailed understanding of the 
beam dynamics in the very complex PSI Ring cyclotron by means of 3D particle simulations. By a 
rough estimation of the initial distribution, according to measurements of beam profile monitors, 
and the time-structure of the beam, realistic simulations of the PSI Ring cyclotron are presented 
and compared to measurements. 

Very good agreement for the radial probe between the simulation and measured data is ob- 
tained by adjusting the injection position, angle, flattop voltage, and the trim coil TC15. These 
parameters are all in agreement with settings obtained from the control room. 

The presented results with a level of detail large enough to predict limiting tails on the extraction 
septum (at beam width levels of 3cr . . . 4(t), and can be seamless extrapolated to future high power 
cyclotrons and enable the precise prediction of crucial parameters, such as losses, based on an 
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existing cw megawatt facility experiences. However a crucial part is the knowledge of the initial 
distribution from the evaluation of measurements of beam profiles and time structure. 

Primary modeling limitations include an accurate knowledge of the initial particle distribution in 
the full 6-D phase space, and the lack of particle-matter interaction in our model. Particle-Matter 
interaction models and resulting struggled primary particles and electrons will play an important 
role when intensity levels increase while at the same time, the losses must be held at present 
levels. We plan to include these effects in future studies, preliminary results on the particle-matter 
interaction model are reported in [15] and ideas for secondary electron creation and field emission 
can be found in [16]. 
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